Gap Opening by Asymmetric Doping in Graphene Bilayers 
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We study the energy gap opening in the electronic spectrum of graphene bilayers caused by asym- 
metric doping. Both substitutional impurities (boron acceptors and nitrogen donors) and adsorbed 
potassium donors are considered. The gap evolution with dopant concentration is compared to the 
situation in which the asymmetry between the layers is induced by an external electric field. The 
effects of adsorbed potassium are similar to that of an electric field, but substitutional impurities 
behave quite differently, showing smaller band gaps and a large sensitivity to disorder and sublattice 
occupation. 
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I. INTRODUCTION 

Graphene materials show unusual structural and 
electronic properties, attracting the attention of the 
scientific community for its potential applications in 
nanoeletronics."^ However, to realize this potential, it is 
important to make graphene semiconducting. Therefore, 
bilayer graphene (BLG), a double layer of graphene with 
AB stacking, is very attractive, since it is a metallic sys- 
tem but it can be made semiconducting by breaking the 
symmetry between the two layer^SEI], More interestingly, 
the gap can be tuned by the magnitude of the symme- 
try breaking, which can be experimentally produced by 
dopin^^ or by an external and perpendicular electric 
fielcP™. 

So far, virtually all theoretical calculations of this effect 
have considered the case of a spatially homogeneous sym- 
metry breaking (with respect to in-plane coordinates), 
which can be experimentally realized by an external elec- 
tric field. However, it is not clear if this approach can 
describe equally well the symmetry breaking induced by 
doping, in which the effects of spatial inhomogeneity, dis- 
order and chemical bonding are present. 

In this work, we address precisely these issues, by a 
combination of ab initio density-functional theory (DFT) 
and tight-binding (TB) calculations. We consider three 
specific situations: boron (acceptor) or nitrogen (donor) 
as substitutional impurities and potassium adsorption 
(donor) in one of the sheets. By comparing our results 
with the case of a homogeneous electric field, we are able 
to understand and to predict the situations in which the 
spatially homogeneous modeling provides a good descrip- 
tion of the gap opening induced by impurity doping. Fur- 
thermore, for the substitutional case, TB calculations in 
large supercells are used to investigate the influence of 
sublattice disorder in the electronic gap. 

Our article is divided as follows: We begin by detailing 
our calculation method in the next Section. Then, in Sec- 
tion |ni| we present our band structure calculations and 



compare the actual doping cases with the homogeneous 
field model for several ordered configurations of impu- 
rities. Next, in Section |IV| we study effects of disorder 
using large supercells within the TB model. Finally, in 
Section |Vj our conclusions are presented. 



II. METHODOLOGY 

Our ab initio calcu lations are based on the density 
functional theorjff^ESl^ pseudopotential method, plane- 
wave basis and supercells with periodic boundary condi- 
tions, as implemented in the the PWSCF package .^^^ We 
work with graphene bilayers with Bernal stacking and 
approximate 10 A of vacuum in the direction perpendic- 
ular to the sheets. The average distance between sheets 
is set to the experimental value for graphite (3.35 A). All 
other geometrical parameters are found by energy min- 
imization until forces are smaller than 10^'^ Ha/Bohr. 
Dipole corrections are employed in order to impose a 
zero electric field in the vacuum region. Hexagonal su- 
percells of different lateral sizes a are used in order to 
vary the dopant concentrations. In this work, we con- 
sider 2x2, 3x3, 4x4, 5x5 and 6x6 supercells. 
Figure [T] shows top views of representative systems for 
a = 3ao (3 x 3) supercells, where oq = 2.46 A is the the- 
oretical lattice constant for BLG. Substitutional (boron 
or nitrogen) and potassium adsorption are shown in Fig. 
[T] (bottom and top panels, respectively). Nitrogen and 
boron substitutions occur in only one of the sheets, thus 
breaking the symmetry between the layers. We choose 
to position the substitutional impurities right on top of 
a carbon atom at the other layer. This choice does not 
affect our main results, as we learn from tests performed 
in other configurations. The Brillouin zone is sampled at 
Nk X Nk X 1 Monkhorst-Pack meshed, in which Nk is 
inversely proportional to the lateral size a, Nk = C/a, 
with C ~ 12ao. For substitutional impurities, a plane 
wave cutoff Ecut = 25Ry is used, along with the gen- 
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eralized gradient approximation (GGA) and Vanderbilt 
ultrasoft pseudopotentials.^^For the K adsorption 
plane wave cutoff of Ecut = 40Ry is used, along the local 
density approximation (LDA) and Von Barth Car (direct 
fit) pseudopotentials with non-linear core correction!^ K 
atoms are placed above the hexagon centers of one of the 
layers. 

For the tight-binding calculations of substitutional im- 
purities, we use the simplest basis set of one Pz orbital 
per carbon atom within the S lonczewsky-Weiss-McClure 
(SWM) parametrization^*'"'^^', with in-plane hopping 70 — 
—2.4 eV, and interlayer hoppings 71 = —0.32 eV {A — A 
coupling) and 73 = —0.30 eV {A — B coupling). The 
onsite matrix elements for carbon atoms is chosen to be 
0.0 eV. These parameters are chosen so as to reproduce 
the DFT band structure of BLG, as shown in Figure [2] 
TB calculations for BLG under a homogeneous electric 
field were performed by adding an energy U to the onsite 
matrix elements of all the atoms in one of the layers, as 
described in detail in the next Section. Finally, a good 
description of the electronic structure of substitutional 
impurities can only be achieved by combining the intro- 
duction of onsite matrix elements at the impurity sites of 
eB — 4.2 eV for boron (the case of nitrogen is completely 
analogous and it was not considered) and a homogeneous- 
field-like U term, as described in Section [rv] 




FIG. 1: Some representative 3x3 supercells studied in this 
work. Upper panel; Potassium adsorption. Lower panel: 
Boron or nitrogen substitutional doping. Black, yellow, light 
blue and red represent carbon (lower sheet), carbon (upper 
sheet), potassium and boron or nitrogen atoms, respectively. 




FIG. 2: DFT (blue full hue) and tight-binding (red dashed 
line) band structures of pristine BLG. 



III. BAND STRUCTURE RESULTS 

Initially, we focus our analysis on the similarities and 
differences between the band structures of BLG in two 
situations: (a) under a homogeneous and perpendicu- 
lar external electric field and (b) under actual doping by 
donor and acceptor species. The DFT band structure 
of BLG under external electric field has been studied by 
Min et al^ and their results show that the simple tight- 
binding model described in the previous Section (includ- 
ing 73 hopping) describes extremely well the evolution of 
the energy gap as a function of the difference in screened 
electrostatic potential between the two layers. Within 
the tight-binding description, this potential difference is 
described by adding an energy U to the onsite matrix 
elements of all carbon atoms in one of the layers, while 
keeping all other parameters fixed. Morever, this param- 
eter U can be easily determined from the band structure, 
since it is precisely equal to the energy gap at the K point 
of the Brillouin zone. Therefore, for simplicity, we adopt 
this tight-binding description for the case of a homoge- 
nous eletric field, but we emphasize once more that it 
provides results that are virtually identical to the DFT 
calculations. 

It is sufficient to analyse in detail the cases of 3 x 3 
and 4x4 supercells, since they are representative of the 
other cases. Figure [3] shows the band structure near the 
Fermi level for different dopants in a 4 x 4 supercell. In 
this case the K point of the BLG Brillouin zone is folded 
into the K point of the supercell. The blue line corre- 
sponds to the actual doping situation while the red line 
is the homogeneous field case, in which U is adjusted to 
reproduce the gap at the K point. One notices that, for 
B and N substitutional doping, the homogeneous field 



FIG. 3: Chemical doping (blue full line) and homogeneous 
electric field (red dashed line) band structures of doped BLG 
in a 4 X 4 supercell. From top to bottom: boron, nitrogen 
and potassium doping. The dashed horizontal line indicates 
the Fermi level. 



FIG. 4: Chemical doping (blue full line) and homogeneous 
electric field (red dashed line) band structures of doped BLG 
in a 3 X 3 supercell. From top to bottom: boron, nitrogen 
and potassium doping. The dashed horizontal line indicates 
the Fermi level. 



description substantially underestimate the band gaps 
and the doping-induced particlc-hole asymmetry is not 
reproduced. The agreement is better for potassium ad- 
sorption, but there are still noticeable differences in the 
band structure. 

Figure [4] shows the band structure near the Fermi level 
for different dopants in a 3 x 3 supercell. Again, the blue 
line corresponds to the actual doping situation while the 
red line is the homogeneous field case. In this particular 
supercell, both the K and K' points of the BLG Bril- 
louin zone are folded into the F point of the supercell, so 
Figure |4] displays twice as many bands as in the 4x4 
case. More importantly, such folding produces degener- 
acy lifting and substantial level repulsion for the case of 
B and N doping. These effects are absent in the homoge- 
neous field approximation, thus making much worse the 
agreement between this approximation and the actual 
doping case. In fact, one notices that gaps are indirect in 
DFT, and this effect is not captured by the homogeneous 
field approximation. Interestingly, the agreement is bet- 
ter for potassium adsorption. The reason is related to 
symmetry: contrary to adsorbed potassium atoms, that 
sit on top of hexagons, B and N substitutional defects 
break the A/B sublattice symmetry, thus coupling more 
strongly K and K points. Therefore, the homogeneous 
field approximation produces a fairly good description of 
the potassium adsorption case. 

For device applications, it is important to describe the 
band gap dependence as a function of defect concentra- 



tion <j — (ag/a)^. This is analogous, in the case of a 
homogeneous electric field, to the gap dependence with 
respect to the field magnitude. The latter case has been 
extensively studied theoretically and it is well known that 
the band gap increases linearly with the field for low field 
values and then saturates at Emax ~ 0.25 eV. In fact, 
from the tight-binding model under the approximation 
of 73 = 0, we expect Eg k, U E,nax / \/ E"^^^ + f/2. As 
an ansatz, we expect that an equivalent or average U 
for the doping case should be proportional to the defect 
concentration, U = aa. Therefore we propose the fol- 
lowing empirical equation to describe the energy gaps as 
a function of defect concentration: 



Figure [5] shows the minimum Kohn-Sham band gap as 
a function of defect concentration for all three dopants. 
Qualitatively, the doped systems behave similarly to case 
of a homogenous electric field and follow roughly Eq. 
[T] showing gap increase and saturation as the concen- 
tration increases. However, there are clear and impor- 
tant differences. First, notice that the gap saturates at 
much lower values (around 0.1 eV) for the substitutional 
dopants boron and nitrogen with respect to potassium 
doping (which in fact follows closely the homogeneous 
field prediction). Second, there are sizeable gap fluctua- 
tions around the empirical curves for the boron and nitro- 
gen cases, probably reflecting the gap sensitivity to the 
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FIG. 5: Band gap as a function of defect concentration for 
boron, nitrogen and potassium doping. Lines are guides to 
the empirical function described in the text. 

particular ordered arrangement of impurities considered: 
most likely these fluctuations are averaged out if a dis- 
ordered arrangement of impurities is considered. These 
effects cannot be captured by a homogeneous field de- 
scription. 

So far, our results show that the homogeneous field ap- 
proximation provides a reasonable description of the gap 
opening in BLG by potassium adsorbates, but it seriously 
fails to describe the behavior of substitutional boron ac- 
ceptors and nitrogen donors. Therefore, it is reasonable 
to question if the gap opening in these substitutional sys- 
tems is robust to a disordered arrangement of impurities, 
in which both A and B sublattices of one of the layers 
are randomly occupied by the defects. We address this 
issue in the following Section. 

IV. DISORDERED ARRANGEMENT OF 
SUBSTITUTIONAL IMPURITIES 

In order to consider a disordered configuration of 
substitutional impurities, one needs to perform calcula- 
tions with very large supercells and configurational aver- 
ages among different realizations. For this purpose, the 
ab-initio approach becomes computationally prohibitive, 
therefore we adopt a tight-binding model. However, our 
tight-binding model must be able to reproduce the DFT 
band structure at least for an ordered configuration of 
substitutional impurities. This is clearly not the case for 
a simple tight-binding model with an homogenous elec- 
tric field, as shown in Section III. 

After an extensive search, we find that a minimum 
tight-binding model that reproduces the DFT band 
structures for ordered arrays of substitutional impurities 



is composed of a short-range and a long-range part. The 
short-range part is described by a change in the onsite 
matrix elements at the impuritity sites. The long-range 
part is modeled by a homogeneous-field-like U onsite po- 
tential for all the atoms belonging to one of the layers 
and it arises from the difference in screened electrostatic 
potential in the two layers caused by the presence of the 
impurities in just one of them. Since our system is com- 
posed of a distribution of screened Coulomb impurities, 
for interdefect distances smaller than the Thomas-Fermi 
screening length in BLCPSl, as it is the case for the range 
of concentrations described in this work, a homogeneous 
field is sufhcient to describe this contribution. It is im- 
portant to note that the long-range contribution is fun- 
damental to open a gap in BLG. We have tested sev- 
eral models containing only short-range modifications in 
the BLG hamiltonian, such as onsite-only or onsite plus 
nearest-neighbor-hopping terms. Such short-range mod- 
els fail to open a gap in BLG. As an example, we show in 
the left panel of Figure |6] the band structure for the 3x3 
structure using two of those models (the pure BLG case 
is also shown for comparison), in which an onsite poten- 
tial of €b = 4.2 eV is used at the impurity site to model 
boron and hopping parameters between the impurity and 
nearest-neighbor carbon atoms are varied between -30% 
and -1-30% from the C-C hopping value, which are rather 
extreme variations. As shown, a gap is not opened. The 
same result is obtained by varying the interlayer hopping 
between boron and carbon and by considering cells with 
different periodicities.-^i- 
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FIG. 6: {left panel) Band structure for pure BLG (black 
full line) and for two short-range models for boron impurities 
in a 3 X 3 cell, using an onsite potential of 4.2 eV at the 
impurity site and nearest-neighbor hoppings of -1.68 eV (red 
dotted line) and -3.12 eV (blue dashed line) between boron 
and carbon, {right panel) DFT (blue full line) and tight- 
binding (red dashed line) band structures of substitutional 
boron-doped BLG in a 3 x 3 supercell. The tight-binding 
model includes long-range and short-range contributions of 
the impurity doping, as described in the text. 

The right panel of Figure |6] shows the DFT and TB 
band structure of an ordered 3x3 array of substitutional 



5 



boron impurities, but the TB model now includes both 
short-range and long-range contributions of the impu- 
rity doping. The DFT (blue) is the same band structure 
shown in Figure |4j The parameters = 4.2 eV and 
U = 0.852 eV are used to fit the TB bands to the DFT 
ones. It is not necessary to vary the hopping elements in 
this shown by the excellent agreement between 

the two curves, as compared with the poor agreement in 
the homogeneous-field model case (Figure |4]). This shows 
that the TB model with long-and short-range contribu- 
tions provides a good description of the effects of substi- 
tutional impurities in graphene and therefore it can be 
used to explore disordered arrangements of such impuri- 
ties. 
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FIG. 7: DOS near the Fermi level for the ordered 3x3 
structure (full red line) and for disordered structures within 
A sublattice (dotted blue line) and within both sublattices 
(dashed black line). DOS for disordered structures involve 
configurational averages for 20 realizations each. In all figures, 
discrete levels are Gaussian-broadened using a 0.02 eV width. 

To achieve this goal, we use supercells with 10,000 sites 
and we distribute 278 boron atoms in one of the layers, 
so as to reproduce the same concentration as the ordered 
3x3 structure. Statistical averaging is performed over 20 
realizations. In order to probe the effects of diferent sub- 
lattice distributions, we consider two different arrange- 
ments: (I) The boron atoms are randomly distributed 
over one of the sublattices only or (2) the boron atoms 
are randomly distributed over both sublattices. The ar- 
rangement (f ) is highly hyphotetical, but it is rather in- 
structive to consider it, as we shall see. Surprisingly, the 
two choices lead to substantially different results for the 
electronic structure, as shown in Figure[7j which displays 
the density-of-states (DOS) near the Fermi level for the 



disordered arrangements (f) and (2), and for the ordered 
3x3 arrangement as well. Notice that the energy gap 
of 0.1 eV for the ordered structure (red full line), con- 
sistent with band structure of Figure [6] If a disordered 
arrangement is imposed to only one of the sublattices 
(blue dotted line) the gap opens up considerably, reach- 
ing almost 0.2 eV. However, and surprisingly, when dis- 
order is imposed to both sublattices (black dashed line) 
the gap shrinks again and shifts to higher energy. This 
situation is more closely related to experimental realiza- 
tions of such doped systems, in which sublattice control 
of impurity positioning is unlikely to be achieved. The 
intrincate dependence of the gap on different types of 
sublattice disorder in the substitutional case in BLG re- 
veals the connection between sublattice symmetry and 
the band structure of graphene systems, as we recall that 
breaking sublattice symmetry is a way to open a gap in 
monolayer graphene. 



V. CONCLUSIONS 

We described the effects of asymmetric doping in the 
band gap of BLG. Adsorbed potassium donors produce 
very similar band gaps, as a function of defect concen- 
tration, as an external and homogeneous electric field 
with varying magnitude. However, the behavior of sub- 
stitutional impurities is quite different, leading to much 
smaller gaps and an acute sensitivity on the disorder and 
sublattice arrangement of impurities. We recall that, in a 
paradigmatic and successful model for donor and accep- 
tor impurities in semiconductors, based on the envelope 
function and effective mass approximations, the short- 
range details of the impurity potential are quite irrele- 
vant. So, one could expect that all donors or acceptors 
(both substitutional and adsorbed) should behave in a 
similar way. Therefore, it is certainly not obvious that 
this is not the case in graphene and we trace the dif- 
ferent behaviors to the local symmetry of both kinds of 
impurities. Therefore, we conclude that models based on 
homogeneous fields can describe quite reasonably the gap 
opening induced by adsorbed donors, but care must be 
taken when the doping is related to substitutional foreign 
species. 
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